/***************************************************************************
 *   Copyright (C) 2006 by Nguyen Quoc Chinh                               *
 *   chinh@pmail.ntu.edu.sg
 ***************************************************************************/
//require two independent variables int l,m for counter

#ifndef MATRIX_MACRO_H
#define MATRIX_MACRO_H
/*==================== for vector */
	//dot product, only for vector[3]
#define DOTPRODUCT(A,B) ((A)[0]*(B)[0]+(A)[1]*(B)[1] + (A)[2]*(B)[2])

	//vector product, only for vector [3]
#define CROSSPRODUCT(c,a,b){(c)[0]=(a)[1]*(b)[2] - ((a)[2]*(b)[1]);(c)[1]=-(a)[0]*(b)[2]+(a)[2]*(b)[0];	(c)[2]=(a)[0]*(b)[1]-(a)[1]*(b)[0];}
	//normalise for vector[3]
#define NORMALISE(A){double sum_0000=sqrt(DOTPRODUCT(A,A)); (A)[0]/=sum_0000;(A)[1]/=sum_0000;(A)[2]/=sum_0000;}

#define DOTPRODUCT_VEC(sum,A,B,dim) {(sum)=0;for(l=0;l<(dim);l++) (sum)+=(A)[l]*(B)[l];}

#define VEC_MUL_VEC(A,B,C,alpha,s0) for(l=0;l<(s0);l++) for(m=0;m<(s0);m++) (A)[l][m]=(alpha)*(B)[l]*(C)[m];

#define VEC_MUL_MAT2(A,B,C,alpha) for(l=0;l<DGR;l++) (A)[l]=(alpha)*DOTPRODUCT((B)[l],(C));

#define VEC_PLUS_VEC(A,B,C,s0,alpha,beta) for(l=0;l<(s0);l++){\
 (A)[l]=(alpha)*(B)[l]+(beta)*(C)[l];}

#define VEC_PLUS_NUM(A,B,s0,num) for(l=0;l<(s0);l++) (A)[l]=(B)[l]+(num);
#define VEC_MUL_NUM(A,B,s0,num) for(l=0;l<(s0);l++) (A)[l]=(B)[l]*(num);
#define VEC_EQU_NUM(A,s0,num) for(l=0;l<(s0);l++) (A)[l]=(num);

/*==================== for 2D matrix */
#define MAT2_MUL_NUM(A,B,num,s0,s1) for(l=0;l<(s0);l++)	for(m=0;m<(s1);m++) (A)[l][m]=(B)[l][m]*(num);
#define MAT2_EQU_MAT2(A,B,s0,s1) for(l=0;l<(s0);l++)	for(m=0;m<(s1);m++) (A)[l][m]=(B)[l][m];

#define MAT2_PLUS_NUM(A,B,num,s0,s1) for(l=0;l<(s0);l++)	for(m=0;m<(s1);m++) (A)[l][m]=(B)[l][m]+(num);
#define MAT2_EQU_NUM(A,num,s0,s1) for(l=0;l<(s0);l++)	for(m=0;m<(s1);m++) (A)[l][m]=(num);
#define MUL_NUM(A,num,s0,s1) for(l=0;l<(s0);l++)	for(m=0;m<(s1);m++) (A)[l][m]*=(num);
#define PLUS_NUM(A,num,s0,s1) for(l=0;l<(s0);l++)	for(m=0;m<(s1);m++) (A)[l][m]+=(num);
#define MAT2_MUL_MAT2(A,B,C,s0,s1,s2) for(l=0;l<(s0);l++) for(n=0;n<(s2);n++){ (A)[l][m]=0;for(m=0;m<(s1);m++) (A)[l][n]+=(B)[l][m]*(C)[m][n];}
#define MAT2_MUL_VEC(A,B,C,s0,s1) for(l=0;l<(s0);l++){(A)[l]=0;for(m=0;m<(s1);m++) (A)[l]+=(B)[l][m]*(C)[m];}
#define MAT2_PLUS_MAT2(A,B,C,s0,s1) for(l=0;l<(s0);l++)	for(m=0;m<(s1);m++) (A)[l][m]=(B)[l][m]+(C)[l][m];
#define PLUS_MAT2(A,B,s0,s1) for(l=0;l<(s0);l++)	for(m=0;m<(s1);m++) (A)[l][m]+=(B)[l][m];

//for converting
#define MAT2_TO_VEC(A,B,s0,s1) for(l=0;l<(s0);l++) for(m=0;m<(s1);m++) (B)[l*(s1)+m]=(A)[l][m];

#define VEC_TO_MAT2(B,A,s0,s1) for(l=0;l<(s0);l++) for(m=0;m<(s1);m++) (A)[l][m]=(B)[l*(s1)+m];

#endif
